# Reproduction files - Study 2 conspiracy myths

# Deepening Bridging and Moving Minds in Stressful Times







# Session Information -----------------------------------------------------




#R version 4.3.2 (2023-10-31 ucrt)
#Platform: x86_64-w64-mingw32/x64 (64-bit)
#Running under: Windows 10 x64 (build 19045)
#
#Matrix products: default
#
#locale:
#  [1] LC_COLLATE=German_Germany.utf8  LC_CTYPE=German_Germany.utf8    LC_MONETARY=German_Germany.utf8 LC_NUMERIC=C                   
#[5] LC_TIME=German_Germany.utf8    
#
#attached base packages:
#  [1] stats     graphics  grDevices utils     datasets  methods   base     
#
#other attached packages:
#  [1] car_3.1-1       carData_3.0-5   patchwork_1.3.0 scales_1.3.0    ggeffects_1.1.5 stargazer_5.2.3 nnet_7.3-18     sjPlot_2.8.12  
#[9] psych_2.2.9     forcats_1.0.0   stringr_1.5.0   dplyr_1.1.4     purrr_1.0.1     readr_2.1.4     tidyr_1.3.0     tibble_3.2.1   
#[17] ggplot2_3.5.1   tidyverse_1.3.2 haven_2.5.3     sjmisc_2.8.9   
#
#loaded via a namespace (and not attached):
#  [1] nlme_3.1-160        fs_1.6.1            lubridate_1.9.2     insight_0.19.0      httr_1.4.4          tools_4.2.2         backports_1.4.1    
#[8] R6_2.5.1            sjlabelled_1.2.0    DBI_1.1.3           colorspace_2.1-0    withr_2.5.0         tidyselect_1.2.0    Exact_3.2          
#[15] mnormt_2.1.1        emmeans_1.8.4-1     compiler_4.2.2      performance_0.10.2  cli_3.6.0           rvest_1.0.3         expm_0.999-7       
#[22] xml2_1.3.3          sandwich_3.0-2      bayestestR_0.13.0   mvtnorm_1.1-3       proxy_0.4-27        minqa_1.2.5         pkgconfig_2.0.3    
#[29] lme4_1.1-31         dbplyr_2.3.0        rlang_1.1.4         readxl_1.4.2        rstudioapi_0.14     farver_2.1.1        generics_0.1.3     
#[36] zoo_1.8-11          jsonlite_1.8.4      googlesheets4_1.0.1 magrittr_2.0.3      Matrix_1.5-1        Rcpp_1.0.10         DescTools_0.99.50  
#[43] munsell_0.5.0       abind_1.4-5         lifecycle_1.0.3     stringi_1.7.12      multcomp_1.4-22     MASS_7.3-58.1       rootSolve_1.8.2.3  
#[50] grid_4.2.2          parallel_4.2.2      crayon_1.5.2        lmom_2.9            lattice_0.20-45     splines_4.2.2       sjstats_0.18.2     
#[57] hms_1.1.2           knitr_1.42          pillar_1.10.1       ggpubr_0.6.0        boot_1.3-28         gld_2.6.6           estimability_1.4.1 
#[64] ggsignif_0.6.4      codetools_0.2-18    reprex_2.0.2        glue_1.6.2          data.table_1.14.6   modelr_0.1.10       nloptr_2.0.3       
#[71] vctrs_0.6.5         tzdb_0.3.0          cellranger_1.1.0    gtable_0.3.1        assertthat_0.2.1    datawizard_0.6.5    xfun_0.37          
#[78] xtable_1.8-4        broom_1.0.3         e1071_1.7-13        coda_0.19-4         rstatix_0.7.2       survival_3.4-0      class_7.3-20       
#[85] googledrive_2.0.0   gargle_1.3.0        timechange_0.2.0    TH.data_1.1-1       ellipsis_0.3.2 




# Calculations ------------------------------------------------------------







# load necessary packages
library(sjmisc)
library(haven)
library(tidyverse)
library(psych)
library(sjPlot)
library(nnet)
library(stargazer)
library(ggeffects)
library(sjPlot)
library(scales)
library(patchwork)
library(car)


# 1. Data preperation, descriptives ---------------------------------------

# read dataset
data_AT <- readRDS("data_at_wght.RDS")


# create new variable: pro health vs. pro liberties
data_AT$pos_pre_group <- car::recode(data_AT$position_pre_cont,
                                     recodes = "1:5 = 0; 6:10 = 1", as.factor = T)

# new variable strong prior opinion
data_AT$strong_prior <- car::recode(data_AT$position_pre_cont,
                                    recodes = "c(1,2,9,10) = 1; 3:8 = 0", as.factor = T)




# 2. Opinion change -------------------------------------------------------

# aggregate opinion change
data_AT <- data_AT %>%
  mutate(abs_position_change = abs(data_AT$position_post_cont - data_AT$position_pre_cont),
         position_change1 = ifelse(abs_position_change <= 0, "no change", "change"),
         position_change2 = ifelse(abs_position_change <=1, "no change", "change"))


# direction of opinion change
data_AT <- data_AT %>%
  mutate(diff_pre_post = position_post_cont - position_pre_cont,
         opinion_formation = ifelse(
           abs_position_change == 0, "Stability", ifelse(
             diff_pre_post > 0 & position_pre_cont>=6, "Polarization", ifelse(
               diff_pre_post < 0 & position_pre_cont>=6, "Depolarization", ifelse(
                 diff_pre_post > 0 & position_pre_cont <= 5, "Depolarization", ifelse(
                   diff_pre_post < 0 & position_pre_cont <= 5, "Polarization", NA))))),
         opinion_formation = as.factor(opinion_formation))



# 4. Regression models ----------------------------------------------------


# Study 2 Model A ---------------------------------------------------------


# set opinion stability as reference category
levels(data_AT$opinion_formation)
data_AT$opinion_formation <- relevel(data_AT$opinion_formation, ref = "Stability")

# levels of treatment group
data_AT$group_treatment <- as.factor(data_AT$group_treatment)
levels(data_AT$group_treatment) <- c("Contestatory", "Collaborative", "Open communication", "Information only", "No information (placebo)")
frq(data_AT$group_treatment)
# reorder factor
data_AT$group_treatment <- factor(data_AT$group_treatment, levels = c("No information (placebo)", "Information only", "Contestatory", "Collaborative", "Open communication"))





# listwise deletion
data_AT <- data_AT %>% 
  filter(!is.na(opinion_formation) & !is.na(group_treatment) & !is.na(strong_prior) &
           !is.na(index_corona_knowledge) & !is.na(leftright) & !is.na(trust_gov) &
           !is.na(trust_exp) & !is.na(index_cognition) & !is.na(index_evaluation) &
           !is.na(index_accuracy))




# Model D -----------------------------------------------------------------

# conspiracy variables
frq(data_AT$conspiracy_myths_1)
frq(data_AT$conspiracy_myths_2) # reversed
frq(data_AT$conspiracy_myths_3)
frq(data_AT$conspiracy_myths_4)
frq(data_AT$conspiracy_myths_5)
frq(data_AT$conspiracy_myths_6)


# rev code conspiracy_myths_2
data_AT$conspiracy_myths_2_rev <- car::recode(data_AT$conspiracy_myths_2, recodes = "1=5; 2=4; 3=3; 4=2; 5=1")



# additive index
data_AT$conspiracy_index <- rowMeans(data_AT[,c("conspiracy_myths_1", "conspiracy_myths_2_rev", "conspiracy_myths_3", "conspiracy_myths_4", "conspiracy_myths_5", "conspiracy_myths_6")], na.rm = T)


# median split > median = hi conspiracy
summary(data_AT$conspiracy_index) # 1.666667

# create grouping variable
data_AT$conspiracy_myths_split <- NA
data_AT[data_AT$conspiracy_index > 1.666667,]$conspiracy_myths_split <- 1
data_AT[data_AT$conspiracy_index <= 1.666667,]$conspiracy_myths_split <- 0



# filter dataset, keep only pre survey conspiracy myth respondents
frq(data_AT$conspiracy_myths_split) # 

data_AT_conspiracy_hi <- data_AT %>% 
  filter(conspiracy_myths_split == 1)

nrow(data_AT_conspiracy_hi) 


# estimate model
at_pol_conspiracy_hi <- multinom(opinion_formation ~ group_treatment + strong_prior + index_corona_knowledge +
                                   leftright + trust_gov + trust_exp + index_cognition +
                                   index_evaluation + index_accuracy, data = data_AT_conspiracy_hi) 

at_pol_conspiracy_hi_N <- nrow(at_pol_conspiracy_hi$residuals) 
at_pol_conspiracy_hi_N # 935 complete cases

# coef plot study 2 model D
coef_plot_study2_mD <- plot_model(at_pol_conspiracy_hi, type = "est", show.values = T, transform = NULL,
                                  auto.label = F, vline.color = "black", 
                                  value.size = 2.5, dot.size = 1.5, line.size = .5, 
                                  value.offset =.35, order.terms = c(1:24),
                                  axis.labels= c("Accuracy motivation", "Neeed for Evaluation",
                                                 "Need for Cognition", "Trust experts",
                                                 "Trust government", "Leftright",
                                                 "Corona knowledge", "Strong prior",
                                                 "Open Communication", "Collaborative",
                                                 "Contestatory","Information Only"))+
  labs(title = "Conspiracy High", caption = paste0("N = ", at_pol_conspiracy_hi_N))+
  theme_sjplot()+
  scale_y_continuous(limits = c(-5,3))
coef_plot_study2_mD[["layers"]][[1]][["aes_params"]]$size <- 0.3
coef_plot_study2_mD[["layers"]][[1]][["aes_params"]]$alpha <- 0.5

# predicted probabilities
pred_at_pol_conspiracy_hi <- ggemmeans(at_pol_conspiracy_hi, terms = "group_treatment", typical = "mean", condition = c(strong_prior =0, gender_male = 0))
pred_at_pol_conspiracy_hi$x <- factor(pred_at_pol_conspiracy_hi$x, levels = c("Open communication", "Collaborative", "Contestatory", "Information only", "No information (placebo)"))

prob_plot_study2_mD <- pred_at_pol_conspiracy_hi %>%
  ggplot(aes(predicted, x))+
  geom_point()+
  geom_segment(aes(x=conf.low, y =x, xend=conf.high, yend=x))+
  facet_grid(~response.level, scales = "free")+
  xlab("")+
  ylab("")+
  theme_bw()+
  scale_x_continuous(n.breaks = 4, labels = label_percent(accuracy = 1))


study2_plotD <- coef_plot_study2_mD / prob_plot_study2_mD






# Model E -----------------------------------------------------------------

# filter dataset, keep only pre survey conspiracy myth respondents
frq(data_AT$conspiracy_myths_split) #

data_AT_conspiracy_lo <- data_AT %>% 
  filter(conspiracy_myths_split == 0)

nrow(data_AT_conspiracy_lo) # 


# estimate model
at_pol_conspiracy_lo <- multinom(opinion_formation ~ group_treatment + strong_prior + index_corona_knowledge +
                                   leftright + trust_gov + trust_exp + index_cognition +
                                   index_evaluation + index_accuracy, data = data_AT_conspiracy_lo) 

at_pol_conspiracy_lo_N <- nrow(at_pol_conspiracy_lo$residuals) 
at_pol_conspiracy_lo_N # 1175 complete cases

# coef plot study 2 model E
coef_plot_study2_mE <- plot_model(at_pol_conspiracy_lo, type = "est", show.values = T, transform = NULL,
                                  auto.label = F, vline.color = "black", 
                                  value.size = 2.5, dot.size = 1.5, line.size = .5, 
                                  value.offset =.35, order.terms = c(1:24),
                                  axis.labels= c("Accuracy motivation", "Neeed for Evaluation",
                                                 "Need for Cognition", "Trust experts",
                                                 "Trust government", "Leftright",
                                                 "Corona knowledge", "Strong prior",
                                                 "Open Communication", "Collaborative",
                                                 "Contestatory","Information Only"))+
  labs(title = "Conspiracy Low", caption = paste0("N = ", at_pol_conspiracy_lo_N))+
  theme_sjplot()+
  scale_y_continuous(limits = c(-3,3))+
  theme(axis.text.y = element_blank())
coef_plot_study2_mE[["layers"]][[1]][["aes_params"]]$size <- 0.3
coef_plot_study2_mE[["layers"]][[1]][["aes_params"]]$alpha <- 0.5

# predicted probabilities
pred_at_pol_conspiracy_lo <- ggemmeans(at_pol_conspiracy_lo, terms = "group_treatment", typical = "mean", condition = c(strong_prior =0, gender_male = 0))
pred_at_pol_conspiracy_lo$x <- factor(pred_at_pol_conspiracy_lo$x, levels = c("Open communication", "Collaborative", "Contestatory", "Information only", "No information (placebo)"))

prob_plot_study2_mE <- pred_at_pol_conspiracy_lo %>%
  ggplot(aes(predicted, x))+
  geom_point()+
  geom_segment(aes(x=conf.low, y =x, xend=conf.high, yend=x))+
  facet_grid(~response.level, scales = "free")+
  xlab("")+
  ylab("")+
  theme_bw()+
  scale_x_continuous(n.breaks = 4, labels = label_percent(accuracy = 1))+
  theme(axis.text.y = element_blank())


study2_plotE <- coef_plot_study2_mE / prob_plot_study2_mE






# Fig. 1 complete
p_part1 <- (coef_plot_study2_mD+coef_plot_study2_mE)


p_part2 <- (prob_plot_study2_mD +prob_plot_study2_mE)
p_part1/p_part2


ggsave("study2_plot_consp_reg_logit.png", height = 7, width = 8, dpi = 300, limitsize = FALSE)




